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Abstract We model the propagation of a coronal shock wave, using nonlinear 
geometrical acoustics. The method is based on the Wentzel-Kramers-Brillouin 
(WKB) approach and takes into account the main properties of nonlinear waves: 
i) dependence of the wave front velocity on the wave amplitude, ii) nonlinear 
dissipation of the wave energy, and Hi) progressive increase in the duration of 
solitary shock waves. We address the method in detail and present results of 
the modeling of the propagation of shock-associated extreme-ultraviolet (EUV) 
waves as well as Moreton waves along the solar surface in the simplest solar 
corona model. The calculations reveal deceleration and lengthening of the waves. 
In contrast, waves considered in the linear approximation keep their length 
unchanged and slightly accelerate. 

Keywords: Magnctohydrodynamics; Waves, Propagation; Waves, Shock 
1. Introduction 

Many solar eruptive events appear to initiate magnetohydrodynamic (MHD) 
waves in the corona. This conjecture is supported by observations of their mani- 
festations. First, these are waves visible in chromospheric spectral lines. For sev- 



eral decades, Moreton waves have been known ( Moreton and Ramsey, 1960 ) ob- 
served in the Ha line. lUchidal (|1968|) proposed that a Moreton wave represented 
the chromospheric trail of a coronal fast-mode wave. More recent studies show 



the possible coronal nature of Moreton waves (see, e.g., Balasub ramaniam, Pevtsov, and NeidigH 
I2007|) . Similar phenomena are propagating waves observed in the He I 10830 
A line (see, e.g., IVrsnak et all I2002bl IGilbert et all l2004j) . Second, type II 
radio bursts are considered to be signatures of coronal shock waves propa- 



gating upwards in the corona (Uchida, 19601. Furthermore, presumable signa 



tures of coronal waves can be observed in soft X-rays (Narukage et al. [2002; 
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microwaves (IWarmuth eFalX [200Ial | White and Thompson! I2005p . and metric 
radio- wavelengths ([Vrsnak et aUl2005p . 

Other phenomena that can be manifestations of MHD waves are large-scale 
wave-like disturbances discovered in 1998 with the Extreme-ultraviolet Imag- 
ing Telescope (EIT) and referred to as "EIT waves" or "extreme-ultraviolet 
(EUV) waves" . This term appears to include phenomena having different phys- 
ical nature and therefore various morphological and dynamical properties. For 
more details, we refer the reader to the Introduction of an accompanying paper 



( Grechnev et al., 2011b I . Note here that describing EUV waves in terms of fast- 
mode MHD waves seems to be possible and correct for those cases when we really 
deal with "wave" phenomena. Just this class of EUV transients is the subject of 
our consideration. 

Thus, fast-mode MHD waves are responsible for a number of solar transients 
and it is important to describe their propagation through the corona. But some 
questions about modeling the waves remain. 

In his original approach, lUchidal (|1968p modeled a coronal disturbance as 
a linear and short fast-mode MHD wave propagating from a point source. To 
calculate a dome-like surface of the wave front and find the position of a Moreton 
wave (as an intersection line of the dome and the solar surface), he used the 
linear geometric acoustics [the Wentzel-Kramers-Brillouin (WKB) approach]. 
Some recent studies (see, e.g., |Wang| 120001 Patsourakos et al.1 I2009[) also used 
the same approach under the assumption of a linear disturbance. In this case, 
the disturbance moves along rays, which curve into regions of a reduced Alfven 
speed. This results in the appearance of wave "imprints" (e.g., Moreton waves) 
running along the spherical solar surface. Note that in the linear approximation 
the amplitude and duration of the wave do not affect the shape of the wave front 
and the speed of its motion along the rays. Neither amplitude nor duration of the 
wave were calculated in the mentioned papers; the geometry of the propagation 
was their only interest. Henceforth in the paper, we consider the wave amplitude 
to be the perturbation magnitude of the plasma velocity in the wave. 

Uchida's model of a linear MHD wave has demonstrated the possibility to 
describe Moreton waves in these terms. Later papers by Uchida with collab- 
orato rs QUchida, Altschuler, and Newkirk 119731 lUchidal 11974ft and recent pa- 



pers ( |Wang| 120001 IPatsourakos et al.[ I2009p . in which the coronal magnetic 
field was calculated from photospheric magnetograms, provided a more accurate 
quantitative description of Moreton and EUV waves in terms of the linear model. 

However, the linear model predicts acceleration of Moreton and EUV waves, 
whereas observations show their systematic deceleration (IWarmuth et a/.ll20011 
2004a). Also, it is often pointed out by many authors that the speeds of coro- 



nal waves sometimes well exceed the fast- mode ones (Narukage et al. [2002; 
IVrsnak et al\ I2002al [Narukage et al\ 12001 IWarmuth et all I2004bl IMuhr et al\ 
I2010p . These facts suggest that the linear approximation does not always cor- 
rectly describe the propagation of those waves. Probably, the disturbance re- 
sponsible for the transients listed above is nonlinear and most likely is a shock 
wave. 

Studies of the propagation of shock waves meet difficulties due to their non- 
linearity. Analytic methods to describe the propagation of shock waves are 
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approximate and often describe the behavior of some extreme classes of nonlinear 
waves (e.g., very strong self-similar waves, or weak waves, etc.). The present 
paper is devoted to the case of a weak shock wave that appears to be the most 
acceptable (see, e.g., IWarmuth et all I2004bl IVrsnak et ail I2002ap . We do not 
discuss the appearance of a shock wave. IVrsnak and Lulicl (I2000p approximately 
described this process, based on an analogy with an acceler ating flat piston. 
The cases of cylindrical and spherical pistons were analyzed by IZic et al\ (|2008p . 
ITemmer et al. \ (|2009p demonstrated how one could describe the kinematics of a 
Moreton wave by using the solution of a simple wave without a discontinuity. In 
our study, we assume that a fast-mode shock wave of moderate intensity appears 
during a solar eruption on the periphery of an active region and decays to a weak 
shock when traveling in the corona. The wave manifests as a Moreton wave and 
an EUV wave on the solar disk. We calculate the propagation of the shock wave 
in terms of the WKB approach taking into account nonlinear effects. 

Such a method in its generally accepted variant involves two independent 
procedures. In the first one, ray trajectories corresponding to the linear approx- 
imation as well as cross sections of ray tubes are calculated. So, the influence 
of the finite wave amplitude on the wave front shape and ray trajectories is 
ignored because in linear acoustics the propagation velocities of disturbances 
are equal to the undisturbed sound speed regardless of their amplitude. In the 
second procedure, the nonlinear variation of the wave amplitude and duration are 
computed along the linear rays obtained. Such a non-self-consistcnt approach is 
fairly useful in some cases, however, nonlinear effects altering the ray pattern and 
the wave velocity disappear completely. Therefore, we develop another method, 
which allows us to consider self-consistently wave propagation and the nonlinear 
variation of wave characteristics. 

The method is described in Section [2] In Section [3] we formulate the problem 
and present results of the analytic propagation model of Moreton and EUV waves 
along the solar surface. Section U contains concluding remarks about the method 
and the results. We note that those who are not interested in the mathematical 
details of the method can read Section 3 without going through Section 2. 



2. Method 

The method of nonlinear geometrical acoustics is based on the method of linear 
geometrical acoustics and allows one to calculate the propagation of disturbances 
with small (but finite) amplitudes through an inhomogeneous medium. The 
linear geometrical acoustics is known to be a method to calculate linear distur- 
bances in the ray approximation (see, e.g., ILandau and Lifshitzl I1987p . In this 
approximation, a solution is found in the form of A (r, t) e 4 *( r '*) where A (r, t) is 
the wave amplitude, and ^ (r, i) is the eikonal, both depending on coordinates 
and time. By substituting this representation for wave perturbations into the 
system of linearized equations of ideal magnetohydrodynamics, one can obtain 
a Hamilton-Jacobi partial differential equation for the eikonal of fast and slow 
magnetosonic waves: 

— + (Vgrad*) + a |grad*| = 0, 
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where V is the undisturbed plasma flow velocity (e.g., the solar wind velocity), 
and a is the magnetosonic speed in plasma. Solving the equation with the method 
of characteristics gives a system of ray equations, which in spherical coordinates 
(r, 9, ip) takes the form QUralova and Uralov, 1994| 

dr k r da 

dt r k dk r ' 

d8 kg da 

dt k dkg ' 



dk r dV r ■ ^ da ^ a ^ 
dt dr r dr kr 



dkg dV r , da , a , 2 „ , dr 

. dk v dV r da ■ m dr , de 

dr a</5 dip dt dt 

where k Tt g^ v are the components of the wave vector k = grad^f, and k is its 
magnitude. The system of equations ((T|) corresponds to the case of a medium 
in steady-state, where only the radial component V r of the undisturbed plasma 
flow exists. By integrating the system ([I}, one can determine the wave front 
shape. 

This approach has been used for modeling coronal fast-mode MHD waves 
dUchida, 1968] [Wang, 2000] ), with only the ray pattern being calculated. How- 



ever, it is essential not only to find the wave geometry, but also to calculate the 
wave intensity. The geometrical acoustics allows the wave amplitude variation 
to be calculated. 

In the linear geometrical acoustics, the energy flux of a disturbance traveling 
in a stationary medium with group velocity qo is directed along the rays, and its 
magnitude is conserved within a ray tube (Blokhintscv, 1981), div(Aeqo) = 0, 



where Ae is the average density of the disturbance energy. In a medium moving 
at a velocity of V, we have to take into account the fact that the wave front 
phase velocity varies as q n = V n + a, with the n index denoting the projection 
normal to the front. The conservation law in this case is div(Aeqg n /a) = 



( Uralov, 1982 |Barnes, 1992[ ) where q = V + qo is the group velocity in a moving 



medium. The average density of the disturbance energy is Ae = p(u 2 + v 2 ) 
where p is the undisturbed plasma density, and u, v are the plasma velocity 
components along the normal to the wave front and across it, respectively. Taking 
into account the relation between the plasma velocity components /1 = v/u 
(Kul ikovsky and Lyubimov, 2005 1, it is possible to relate the variation of the 



wave amplitude to the normal cross section dS of the ray tube formed by a 
bundle of rays: 

dSqpu 2 (l + /J 2 ) — = const. (2) 
a 

Thus, the general approach to determine the wave amplitude in the ray method 
relies on calculating the cross section of the ray tube. 
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There are various techniques for calculating cross sections. They are discussed 
in lKravtsov and Orlovl(|1990p . We calculate cross sections by using the Jacobian 
of the transformation to ray coordinates. The volume element dW of the ray 
tube is expressed in terms of ray coordinates (771 ,r]2,t) as 

dW = dxdydz = r 2 sin 9drd9dip = r 2 s'mOD (t) drjid^dt, 

where D(t) is the Jacobian of the transformation from spherical coordinates to 
ray ones. Then for the cross section of the ray tube, dS, we have 



, C dW 2 • 

db = — — = r sinfl drjidrj 



2- 



da q 
where a is the ray tube length. Substituting it into ^ yields 

D (t) r 2 sin 6pu 2 (l + u 2 ) — = const. (3) 



To calculate the Jacobian, we use a method based on numerical integration 
of a so-called adjoint system of equations. This system consists of differential 
equations for the derivatives dr/drji^, dO/drji^, dtp/drji t 2, dkr/drji^, dkg/drji^, 
dk v / 1 drj\ t 2 and is derived from ([1} by differentiating the equations with respect 
to ray coordinates r\\ and 772. Let the ray coordinates 771 and 772 be the angles 
defining the direction of the outgoing ray from a point source at the initial 
moment. Note that the ray coordinates do not need to be explicitly determined 
when the adjoint system is being derived. This becomes essential to specify 
the initial values for the desired functions. For the case considered, the adjoint 
system has the following form (in view of the symmetry about 771 and 772, instead 
of 12 equations we give only six for one variable 77): 

d f 9r\ dV r k r da a dk r ak r dk da dk h®( ® a 
dt \dn ) dn k dn k dn k 2 dn dk r dn dn \ dk r 



dr d6 d [ d6\ kg da a dkg akg dk da dk , d ( d> 



r- 



+ tt; ^ + k- 



dn dt dt \dn J k drj k dn k 2 dn dkg dn dn \ dk 

dr . dip n d8 dip , d ( dip 

— sine— +rcos9—— + rsm9— — 
dn dt dn dt dt \ dn 

dak, ^ a dk, ak v dk da dk ^^d / da 

dn k k dn k 2 drj dk, drj dr] \dk, 

d ( dk r \ dV r dk r d ( dV r \ dadk d ( da\ dakg + k 2 

k 



dt \ dn J dr drj drj \ dr J dr drj drj \dr J drj kr 
a 



< dk a(k 2 + kl) dr | a / gfa | ^ dk, 

k 2 r v drj kr 2 dn kr \ drj v drj 

dr dkg d ( dkg \ dV r dk r d ( dV r \ dadk d f da 

d^j~dt +r di \~drj) = ^~d9~d^ ~ r £hj \~d9 J ~ d9~(ki dy\d6 

dak 2 cot 9 ak 2 dk 2ak, dk, 

+ a 7 FT cot "I 1 cot 0— ^ - 

dn k k z dn k drj 
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gfc; 89 dkgdr d fdr\ 

ksm 2 9dr] dr] dt 6 dt \8r)J ,[ ' 

dr . dk v n d9 dk v . n d ( dk v 

— smd—rf- + r cos 9— —rf- + rsm9— — ^ 
or] dt or] dt dt \ or] 

_9K9Av _ k d_ fdVA _ dadk_ k _d_ f da\ _ Q , dO_ k dr 

dp dr] r dr] \ dip J dp dr] drj \dp J dr] v dt 



, dk w dr . „ , d ( dr\ dk w d9 , dr ,d9 

r cos k v — cos — 

dt dr] dt 

86 d9 , a d (89 
+k v rsm0— — - k v r cos 0— — 
or] dt dt \ or] 



where 



8V r 8V r dr a da ^—^ ( da dr a da dk a \ dk 1 ^ <9fc Q 

9?/ ^— ^ <9r Q (9?; ' 9?/ ^— ^ \9r a 9?/ dk a dr] ) ' <9?7 ^— ^ a 9?/ ' 

d ( da \ ^ Z' ^ ( ® a \ ^ ra ® ( ® a \ dk(% 
dr] \dk/3 J ^ \dr a \dkp J dr] dk a \dkp J dr/ 

d ( da\ >— ^ / d ( da\ dr a d ( da \ dk a 
dr] \drp J ^ \dr a \drp J dr] dk a \drp J dr] 

^(^)=E(^(|^)^)' r a ,e = {r,9,<p}, k a , = {k r ,ke,K}. 

Thus, to calculate the propagation of a linear wave and its intensity, at first we 
have to integrate numerically the ray equations system (|T|) and adjoint system (0]) 
and then determine the amplitude variations by means of J3]). 

A nonlinear flat disturbance in an ideal homogeneous medium is described by 
a simple wave solution and propagates at the supersonic speed determined by its 
amplitude ( |Kulikovsky and Lyubimov, 2005| ) . A fast-mode simple- wave element 
with plasma velocity component u normal to the front moves at a + ku, where 
K = (1/a) (d (pa) I dp) is the numerical coefficient depending both on the plasma 
beta and the angle between the wave vector and the magnetic field. We do not 
give here the bulky explicit expression for n. We only note that values of k are 
restricted by the limits: 1/2 < n < 3/2. The fact that each simple-wave element 
travels at its own speed causes wave profile deformation and the appearance 
of a discontinuity. If a moderate amplitude simple wave has a triangular profile 
before the discontinuity appears, it will take the shape of the right-angled triangle 
after the discontinuity forms, with the discontinuity being the leading edge of 
the disturbance. Note that any nonlinear disturbance profile of a finite duration 
tends asymptotically to this shape. Let U s h be the jump of the plasma velocity 
component u in the discontinuity. Then, in the nonlinear geometrical acoustics 
approximation, the discontinuity moves at a speed of a + nU s h/2- Taking into 
account this increase in the wave front speed in the ray equations, we are able to 
correctly describe the propagation of weak shock waves. Then the ray equation 



SOLA: modeling_waves.tex; 19 January 2013; 23:36; p. 6 



Modeling MHD shock wave propagation 



system ([T]) becomes ( jUralova and Uralov, 1994| 

dr f KU sh \ k r da 

dt = Vr+ { a + —)T + k dk r ' 

dO ( KU sh \ k e da 

T Tt = \ a + — )T + k dh- 

dk r dV r da a , 2 2 , 

~dt = ~~bT kr ~dr' k+ kr'( ke + k ^> 

dkg dV r da a 2 dr 

. dk v dV r da ■ m dr , Q de 

rsmO—- 1 - = — - — k r — ——k — smOku,— k^rcosu—. 

dt dip d(p dt dt 

The generally accepted method of nonlinear geometrical acoustics ignores the 
additional term KU s h/2 as it is small. However, it is the term that is responsible 
for wave deceleration due to the amplitude damping. Besides, estimations within 
the framework of the perturbation theory suggest that the ray pattern variation 
due to nonlinearity is a correction of the same order of magnitude as the non- 
linear variation of the wave amplitude is. It is therefore important to take this 
into account in the nonlinear geometrical acoustics approximation. 

Ray equation system is not closed now because it includes the wave 
amplitude. In the linear approximation, an amplitude variation can be deter- 
mined from ((3]). The nonlinear wave amplitude undergoes additional damping 
associated with energy dissipation in the discontinuity. As the amplitude, we 
take the value of the jump U s h- Variations of amplitude U s h and duration T s h 
of a weak shock wave having a triangular compression phase may be calculated 
as ( jUralov, 1982] ) 

T l \ rn mli T l \ " T l KU 1 



u " = *{ 1 + £) ' T ° h = T *{ 1 + i) ' dT = -£f (6) 

where t\ is the duration increment of the simple wave with an amplitude of u\ ; 
T* is the initial duration of the disturbance. Note that the laws (|6]) of the weak 
shock wave damping are derived by using values of the amplitude and duration 
of a simple wave, from which the discontinuity forms. The value of u\ can be 
determined from the expression similar to . 

Thus, solving numerically of 19 ordinary differential equations (|5|), and © 
enables us to compute propagation of a weak shock wave in an inhomogencous 
medium, its amplitude and duration. 



3. Analytical Modeling of Wave Propagation 

In this section, we employ the nonlinear geometrical acoustics method to de- 
scribe the propagation of large-scale wave-like transients, namely EUV and 
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Figure 1. The Alfven speed distribution in the solar corona model 10, JS). 



Moreton waves. With respect to "EUV wave" phenomena, we address only those 
disturbances that are associated with a fast-mode MHD shock wave. 

To calculate the propagation of a shock wave, we have to specify the solar 
corona model as well as position and parameters of the shock wave at the 
initial moment. We use a simple hydrostatic corona model to demonstrate the 
main particularities of the method and compare results with those obtained 
in the linear approximation. The corona is considered to be isothermal with 
temperature T — 1.5 x 10 6 K and sound speed c = 181 km s _1 . The plasma 
density is distributed in accordance with the barometric law (for details see, e.g., 
IMann et oil [TOM]) : 

n(r)=T*«p(§(^-l)), (7) 

where no = n(i? Q ) = 3 x 10 s cm~ 3 is the density at the base of the corona, Rq 
the solar radius, H = 2R gas T /rhMngo ~ 70 Mm the density scale height, g® 
the acceleration of gravity on the solar surface, Mh the molar mass of hydrogen, 
fh = 1.27 the average atomic weight of an ion, and R gas the gas constant. Let us 
assume that we have a magnetic field with only a radial component (for details 
see, e.g.. IMann eTal! 12003"]) : 

B r = ±B (^j , (8) 

where Bq = 2.3 G is its value at the base of the corona. The sign in © depends 
on the solar hemisphere, but it is not meaningful here. The same model was 
applied by lUchidal (fT0rJ8) . 

In the corona model ([7]) and ©, the Alfven speed increases with height, 
peaking at Rq/AH = 2 A3 Rq (Figure [I]). Refraction makes ray trajectories 
curved towards regions of the lower Alfven speed. The solar corona has therefore 
waveguide properties. A portion of the wave energy flux is captured by the 
coronal waveguide and propagates along the solar surface, giving rise to an EUV 
and Moreton wave. In the treatment used, the EUV front is observable due to 
the plasma compression produced by the coronal shock wave. Since the plasma 
density decreases rapidly with height, a plasma layer near the solar surface 
contributes substantially to the EUV front emission. The layer thickness is about 
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Figure 2. Propagation of a coronal shock wave as seen in a 3D image (left) and a 2D cross 
section with ray trajectories (right). The solid arcs drawn across the rays represent shock wave 
fronts. The dashed line marks a height of 60 Mm above the solar surface. The wave front 
velocity along this line corresponds to the velocity of an EUV wave. The wave source height 
is 80 Mm. 



the density scale height H . So, to estimate the EUV front position, we have to 
find the intersection line of the calculated shock front and a spherical surface of 
radius < (Rq+H). The Moreton wave corresponds to the chromospheric trail of 
the coronal shock wave, i.e. it can be found at the intersection line of the shock 
front and the upper chromosphere. 

For modeling, it is important to assign the initial characteristics of the shock 
wave and its start position. We have to specify the initial duration (or length) 
of the wave and its initial amplitude on some surface. In this study, these values 
are given according to the strong point-like explosion theory. The wave source 
located at a height of 80 Mm is characterized by the energy e whose release 
produces the shock wave. When the wave covers a distance of A = (e/ p*^) 1 / 3 , 
with and a* being respectively the plasma density and the fast-mode speed 
at the explosion point, the compression phase profile of the wave is assumed to 
be triangular. The compression phase length is equal to A and the amplitude is 
X a *i with x being a coefficient of the order of unity. In this paper, we employ a 
value x = 1 except for the calculations given in Figure [7J The initial duration 
of the compression phase is supposed to be equal to T» = A(k)/a*(k). We 
believe that a shock wave arises on the periphery of an active region located 
within the explosion cavity A. This manner to assign initial values does not rely 
on the specific mechanism of the wave initiation. It is essential only that the 
energy release producing a shock wave is impulsive. For instance, a wave can be 
produced by a compact piston acting for a short term (an abruptly accelerating 
filament and its magnetic envelope). 

Figures|5H3present the results of our modeling. Figure[5]illustrates a 3D image 
of the shock wave front and the respective 2D section including ray trajectories. 
The rays go out from the initial surface of size A. The wave front gets inclined 
over the solar surface, with its inclination increasing in time. 

Figure [3] shows the distance-time plots of EUV waves (a) as well as the time 
plots of their velocities (b) along the surface 60 Mm (relative to the solar surface). 
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Figure 3. The distance-time plot of EUV waves (a) and the time evolution of their velocities 
(b) along the surface at 60 Mm height for different energies e specified in panel b. The lowest 
(dotted) line represents the linear EUV wave. 
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Figure 4. The kinematic plots of a Moreton wave (thick) and an EUV wave (thin), which are 
produced by a single coronal wave in the linear (dotted) and nonlinear (solid) consideration. 

The curves are given for different values of the wave source energy e or different 
values of the initial shock wave length, as follows from A = (e/p^aj) 1 / 3 . The 
EUV wave velocity decreases appreciably due to the nonlinear damping of the 
coronal wave amplitude. After the amplitude has substantially decreased, the 
shock wave propagates as a linear one. Having reached its minimum, the wave 
velocity slightly increases. This is due to the shock front inclination over the 
solar surface that becomes more and more significant with time and associated 
with the waveguide properties of the lower quiet Sun's corona. The larger the 
wave front inclination, the higher the Moreton and EUV wave velocity. Note 
that calculated Moreton and EUV wave acceleration must be difficult to observe 
since it occurs when the wave amplitude becomes low (see Figure [5]). 

Figure |4] shows positions (a) counted along the surface and velocities (b) of a 
Moreton wave (thick lines) and an EUV wave (thin lines), which are produced 
by a single coronal wave. The dotted lines correspond to a linear coronal wave 
and the solid lines are for a shock one. The linear Moreton wave velocity is lower 
because the Alfven speed at the corona base is smaller than that at a height of 
60 Mm. Note that a Moreton wave decelerates even in the linear case. This fact 
is associated with the initial wave source position at a height of 80 Mm above 
the photosphere. The EUV wave in the linear case does not decelerate since the 
wave source is located roughly at the same height that the EUV wave is (a little 
bit higher the dashed line in Figured]). 
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The increase in the propagation velocity of the wave front with height results 
in a front inclination relative to the solar surface. So one can observe the wave 
signatures corresponding to the different heights {e.g., EUV waves and Morcton 
waves) to be shifted. Increase of the front inclination with time determines a 
value of this shift. Besides, its time evolution is also determined by the observer 
position and the sight angle. Such a consideration demonstrates the possibility 
to explain the offset between Moreton and EUV waves as well as the Hel-Ha 
offset observed bv lVrsnak et aZ.I (|2002b|) since waves in the He I 10830 A line are 
similar morphologically to EUV ones. 

Figure [5] presents time dependence of the EUV wave amplitude for different 
values of the wave source energy e. A disturbance having higher e and greater 
length decays more slowly as follows from ([6]) and as is seen in Figure [5j The 
Moreton wave amplitude varies in a similar manner, but it has smaller values. 

Another effect associated with the nonlinearity of an EUV wave is the increase 
in its duration T s h and, respectively, length L = aT s h (also referred to as the 
wave profile broadening). In the linear case, the wave duration is constant (under 
the assumption of a steady-state medium). Figure |6] gives the time evolution of 
the ratio of the EUV wave length to its initial size A. If initial amplitudes of 
disturbances are equal, the relative extension will be faster for disturbances 
having the lower initial energy (and the shorter initial length) . 

With respect to the damping of shock waves having the same initial lengths, 
amplitude decrease is faster for a wave with a higher amplitude. Therefore, 
shock waves with different initial amplitudes decay to the same level after ap- 
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Figure 7. The amplitude damping of shock waves (solid lines) and linear ones (dashed lines) 
calculated for the wave source energy e = 10 29 erg and different initial amplitudes. The initial 
amplitude values are specified in the figure as ratios x °f those to a fast-mode speed in the 
wave source. 



proximately equal duration (FigureEl solid curves). For comparison, we also plot 
the amplitude curves of linear disturbances, which hold an initial ratio of the 
wave amplitudes throughout propagation. The mentioned property of nonlinear 
waves (in Figure [7]) allows us to be not precise about the value of the initial 
wave amplitude. Therefore, we can use the value that follows from the strong 
point-like explosion theory. 



4. Discussion and Conclusion 

We have modeled the propagation of shock-associated EUV waves and Moreton 
waves, using the nonlinear geometrical acoustics method. This method takes into 
account characteristic properties of nonlinear waves: i) dependence of the wave 
velocity on its amplitude, ii) wave energy dissipation in the shock front, and 
Hi) wave duration increase with time. The method allows one to calculate the 
nonlinear evolution of the shock wave and its propagation pattern. However, the 
generally accepted variant of this approach includes nonlinearity effects only for 
describing the amplitude and duration of a wave, but not for the wave velocity 
value. So, using such a non-self-consistent approach results in the loss of an 
important effect concerning wave kinematics. 

We have applied another approach, appending an additional term to the ray 
equations. This has allowed the finite wave amplitude to be taken into account. 
We have solved self-consistently the modified ray equations and the equations 
describing the wave amplitude and duration evolution along a ray tube. It is this 
approach that has been developed in this paper to analyze coronal shock wave 
propagation along the solar surface. 

One of the results of our analysis is the deceleration of EUV and Moreton 
waves at the initial stage of propagation. Since we use a spherically symmetric 
and isothermal model of the solar corona, deceleration is a direct consequence of 
their nonlinearity. Thus, EUV and Moreton waves having a sufficient amplitude 
(and therefore being observable) have to decelerate in the quiet Sun's regions 
where average plasma parameters are constant along the solar surface. Note that 
the large-scale waves under study are registered just in these regions. Calculated 
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wave deceleration is supported by the EUV and Moreton wave observations 
analyzed bv lWarmuth et aO(|2001L[2004ap . Also we notice here that EUV waves 
not associated with a coronal MHD wave show only slight or no deceleration 



(see, e.g.,|Wills-Davey and AttHIH \M)§$ . 



The simple corona model also let us find other features of wave kinematics, 
e.g., i) the wave source height has effect on the initial portion of the velocity 
plot, and ii) the rate of wave deceleration and damping becomes lower as the 
wave source energy (or the wave length, respectively) grows. All these findings 
are applied for modeling EUV wave propagation in the 17 January 2010 event 
in an accompanying paper ( jGrechnev et al., 201 lap . 

Modeling waves in the linear approximation does not reveal their deceleration. 
On the contrary, linear waves undergo only small acceleration caused by a slightly 
increasing inclination of the coronal wave front over the solar surface. This effect 
was first discovered bv lUch ida (1968;) . However, because of the error that he made 
in the expression for the barometric distribution of the coronal plasma density 
(the scale height was halved), the wave front inclination over the solar surface 
was very large. As a result, linear waves underwent considerable acceleration in 
Uchida's original model. 

Another important result of our modeling is the duration (and length) increase 
of Moreton and EUV waves. This effect is also confirmed by observations (see, 
e.q., IWarmuth et al. 11200 ll|Veronig et aL|l2010p . In contrast, a linear disturbance 
keeps its duration unchanged in a steady-state medium. Note that in the linear 
approximation, the wave amplitude and duration also vary due to the viscosity, 
the thermal conductivity and the finite plasma conductivity, however, these 
effects are negligible against nonlinear factors. 

To summarize, we believe that wave deceleration and its duration increase, 
both being the attributes of shock wave evolution, point out the crucial role of 
nonlinearity in the behavior of EUV and Moreton waves (at least, it concerns 
some of them). 

In conclusion, we will briefly discuss the method limitations for solving the 
shock wave propagation problem. The main limitation is associated with laws Q 
of a shock wave damping, which are derived by using the relations for simple flat 
MHD waves in a homogeneous medium. So, we have to meet two requirements. 
First, the shock wave length should be smaller than the radius of curvature of the 
wave front and the smallest medium variation scale. The fulfilment of these con- 
ditions also ensures validity of the linear ray approximation ((1|), which involves 
actually even less limitations. The smallest variation scale in our modeling is that 
of plasma density. So, we are aware that our computation lies at the boundary of 
applicability of nonlinear geometrical acoustics since characteristic shock wave 
length A and density scale pj |Vp| arc of the same order of magnitude. 

Second, the nonlinear factor U s h/a should be small. Under this condition, 
damping laws ^ are derived and this very condition ensures a correct cal- 
culation of the terms involving U s h in the ray equations system ([5]). With 
respect to the limitations of laws (O, the point-like atmospheric explosion theory 
dKestenboim, Roslyakov, and Chudov, 1974 ) suggests that these laws describe 
satisfactorily spherical shock wave propagation up to U s h/a < 1. When we choose 
the initial value U s h — a in Section [31 we do not therefore go beyond the scope 
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of the application of relations ©. However, ray equations ([5]) and equations (0J 
at U s h/a ps 1 are able to yield an error in calculations. Nevertheless, this error 
is insignificant due to the nearly spherical shape of the wave front at the initial 
phase of propagation and then it disappears owing to rapid decrease in U s h ■ 
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